Lists

rm(list=ls())

Defining the data I/O directory for the Eddy Master File

**

path_eddy<-"C:\\Users\\Tommy\\flux\\Data-exploring\\02_Concord\\"
path.in_eddy<-paste(path_eddy,"01_Proccessed_Data",sep="")
path.out_eddy<-paste(path_eddy,"03_combined_data",sep="")
ver<-"Master_Eddy" 
file.name_eddy<-paste("master_eddy_pro_concord",sep="")

read in eddypro full output Master file, parse variable names

data_master_eddy<-read.csv(paste(path.in_eddy,"\\",ver,"\\",file.name_eddy,".csv",sep=""),
                  header=F,
                  skip=3,
                  na.strings=c(-9999),
                  stringsAsFactors = F)
colnames(data_master_eddy)<-colnames(
  read.csv(paste(path.in_eddy,"\\",ver,"\\",file.name_eddy,".csv",sep=""),
           header=T,
           skip=1))

data_master_eddy

Defining the data I/O directory for the Master met_data

path_met<-"C:\\Users\\Tommy\\flux\\Data-exploring\\02_Concord\\"
path.in_met<-paste(path_met,"01_Proccessed_Data",sep="")
path.out_met<-paste(path_met,"03_combined_data",sep="")
ver<-"met_data" 
file.name<-paste("MET_data_master",sep="")

#read in Met_Data Master file, parse variable names and define N/As#

met_data_master<-read.csv(paste(path.in_met,"\\",ver,"\\",file.name,".csv",sep=""),
                  header=F,
                  skip=4,
                  na.strings=c("NAN"),
                  stringsAsFactors = F)
colnames(met_data_master)<-colnames(
  read.csv(paste(path.in_met,"\\",ver,"\\",file.name,".csv",sep=""),
           header=T,
           skip= 1))

met_data_master

NA

Parsing the time stamp converting it into a POSIXlt vector

interpreting date and time into new timestamp column

Then taking that time stamp column and turning each time into a unique number (time.id) so I can join based on that. As it can be really tricky to join/merge based on time stamps alone

Or I could make sure both time stamps are characters and match them that way

Finally ploting time.id to make sure my times translate linearily

data_master_eddy$TIMESTAMP<-strptime(paste(data_master_eddy$date,data_master_eddy$time,sep=" "),format="%m/%d/%Y %H:%M", tz = "GMT")

data_master_eddy$time.id<-data_master_eddy$TIMESTAMP$year+1900+(data_master_eddy$TIMESTAMP$yday)/366+(data_master_eddy$TIMESTAMP$hour)/366/24+ (data_master_eddy$TIMESTAMP$min)/366/24/60

data_master_eddy$time.id[1:50]

 [1] 2019.450 2019.450 2019.450 2019.450 2019.450 2019.450 2019.450 2019.450 2019.451
[10] 2019.451 2019.451 2019.451 2019.451 2019.451 2019.451 2019.451 2019.451 2019.451
[19] 2019.451 2019.451 2019.451 2019.451 2019.451 2019.451 2019.451 2019.452 2019.452
[28] 2019.452 2019.452 2019.452 2019.452 2019.452 2019.452 2019.452 2019.452 2019.452
[37] 2019.452 2019.452 2019.452 2019.452 2019.452 2019.452 2019.452 2019.453 2019.453
[46] 2019.453 2019.453 2019.453 2019.453 2019.453

plot(data_master_eddy$time.id)

which(duplicated(data_master_eddy$time.id))

integer(0)

#Taking the met_data and turning the time stamp into posixt format#

met_data_master$TIMESTAMP<-strptime(met_data_master$TIMESTAMP,
                           format ="%m/%d/%Y %H:%M", tz = "GMT")

met_data_master$TIMESTAMP[1:20]

 [1] "2019-06-25 09:00:00 GMT" "2019-06-25 09:30:00 GMT" "2019-06-25 10:00:00 GMT"
 [4] "2019-06-25 10:30:00 GMT" "2019-06-25 11:00:00 GMT" "2019-06-25 11:30:00 GMT"
 [7] "2019-06-25 12:00:00 GMT" "2019-06-25 12:30:00 GMT" "2019-06-25 13:00:00 GMT"
[10] "2019-06-25 13:30:00 GMT" "2019-06-25 14:00:00 GMT" "2019-06-25 14:30:00 GMT"
[13] "2019-06-25 15:00:00 GMT" "2019-06-25 15:30:00 GMT" "2019-06-25 16:00:00 GMT"
[16] "2019-06-25 16:30:00 GMT" "2019-06-25 17:00:00 GMT" "2019-06-25 17:30:00 GMT"
[19] "2019-06-25 18:00:00 GMT" "2019-06-25 18:30:00 GMT"

#Making sure timestamp columns line up#

met_data_master$TIMESTAMP[1:10]

 [1] "2019-06-25 09:00:00 GMT" "2019-06-25 09:30:00 GMT" "2019-06-25 10:00:00 GMT"
 [4] "2019-06-25 10:30:00 GMT" "2019-06-25 11:00:00 GMT" "2019-06-25 11:30:00 GMT"
 [7] "2019-06-25 12:00:00 GMT" "2019-06-25 12:30:00 GMT" "2019-06-25 13:00:00 GMT"
[10] "2019-06-25 13:30:00 GMT"

data_master_eddy$TIMESTAMP[1:10]

 [1] "2019-06-14 17:30:00 GMT" "2019-06-14 18:00:00 GMT" "2019-06-14 18:30:00 GMT"
 [4] "2019-06-14 19:00:00 GMT" "2019-06-14 19:30:00 GMT" "2019-06-14 20:00:00 GMT"
 [7] "2019-06-14 20:30:00 GMT" "2019-06-14 21:00:00 GMT" "2019-06-14 21:30:00 GMT"
[10] "2019-06-14 22:00:00 GMT"

met_data_master

#creating a time id for the MET Data so I I can join the MET and Eddy Pro Data#

met_data_master$time.id <-met_data_master$TIMESTAMP$year+1900+(met_data_master$TIMESTAMP$yday)/366+(met_data_master$TIMESTAMP$hour)/366/24 + (met_data_master$TIMESTAMP$min)/366/24/60 

met_data_master$time.id[1:20]

 [1] 2019.479 2019.479 2019.479 2019.479 2019.479 2019.479 2019.480 2019.480 2019.480
[10] 2019.480 2019.480 2019.480 2019.480 2019.480 2019.480 2019.480 2019.480 2019.480
[19] 2019.480 2019.480

plot(met_data_master$time.id)

which(duplicated(met_data_master$time.id))

integer(0)

met_data_master

#Joining the Met_Data and Eddy Pro Data Sets#

with merge

combo_master_ed_met<- merge(met_data_master[,-which(colnames(met_data_master)=="TIMESTAMP")], data_master_eddy[,-which(colnames(data_master_eddy)=="TIMESTAMP")], by = "time.id")

combo_master_ed_met

colnames(combo_master_ed_met)

  [1] "time.id"                       "RECORD"                       
  [3] "BattV_Avg"                     "PTemp_C_Avg"                  
  [5] "AM25T_ref_Avg"                 "TC_Avg.1."                    
  [7] "TC_Avg.2."                     "TC_Avg.3."                    
  [9] "TC_Avg.4."                     "TC_Avg.5."                    
 [11] "TC_Avg.6."                     "TC_Avg.7."                    
 [13] "TC_Avg.8."                     "TC_Avg.9."                    
 [15] "TC_Avg.10."                    "TC_Avg.11."                   
 [17] "TC_Avg.12."                    "TC_Avg.13."                   
 [19] "TC_Avg.14."                    "TC_Avg.15."                   
 [21] "TC_Avg.16."                    "TC_Avg.17."                   
 [23] "TC_Avg.18."                    "TC_Avg.19."                   
 [25] "TC_Avg.20."                    "TC_Avg.21."                   
 [27] "TC_Avg.22."                    "TC_Avg.23."                   
 [29] "TC_Avg.24."                    "TC_Avg.25."                   
 [31] "AirT_Avg"                      "RH_Avg"                       
 [33] "AtmPressure_Avg"               "NR_mV_Avg"                    
 [35] "NR_Wm2_Avg"                    "PAR_in_mV_Avg"                
 [37] "PAR_in_uEm2_Avg"               "PAR_out_mV_Avg"               
 [39] "PAR_out_uEm2_Avg"              "SHF_1_mV_Avg"                 
 [41] "SHF_1_Wm2_Avg"                 "SHF_2_mV_Avg"                 
 [43] "SHF_2_Wm2_Avg"                 "WaterP_Avg"                   
 [45] "WaterT_Avg"                    "Precip_mm_Tot"                
 [47] "VWC_Avg"                       "EC_Avg"                       
 [49] "T_Avg"                         "P_Avg"                        
 [51] "PA_Avg"                        "VR_Avg"                       
 [53] "filename"                      "date"                         
 [55] "time"                          "DOY"                          
 [57] "daytime"                       "file_records"                 
 [59] "used_records"                  "Tau"                          
 [61] "qc_Tau"                        "H"                            
 [63] "qc_H"                          "LE"                           
 [65] "qc_LE"                         "co2_flux"                     
 [67] "qc_co2_flux"                   "h2o_flux"                     
 [69] "qc_h2o_flux"                   "H_strg"                       
 [71] "LE_strg"                       "co2_strg"                     
 [73] "h2o_strg"                      "co2_v.adv"                    
 [75] "h2o_v.adv"                     "co2_molar_density"            
 [77] "co2_mole_fraction"             "co2_mixing_ratio"             
 [79] "co2_time_lag"                  "co2_def_timelag"              
 [81] "h2o_molar_density"             "h2o_mole_fraction"            
 [83] "h2o_mixing_ratio"              "h2o_time_lag"                 
 [85] "h2o_def_timelag"               "sonic_temperature"            
 [87] "air_temperature"               "air_pressure"                 
 [89] "air_density"                   "air_heat_capacity"            
 [91] "air_molar_volume"              "ET"                           
 [93] "water_vapor_density"           "e"                            
 [95] "es"                            "specific_humidity"            
 [97] "RH"                            "VPD"                          
 [99] "Tdew"                          "u_unrot"                      
[101] "v_unrot"                       "w_unrot"                      
[103] "u_rot"                         "v_rot"                        
[105] "w_rot"                         "wind_speed"                   
[107] "max_wind_speed"                "wind_dir"                     
[109] "yaw"                           "pitch"                        
[111] "roll"                          "u."                           
[113] "TKE"                           "L"                            
[115] "X.z.d..L"                      "bowen_ratio"                  
[117] "T."                            "model"                        
[119] "x_peak"                        "x_offset"                     
[121] "x_10."                         "x_30."                        
[123] "x_50."                         "x_70."                        
[125] "x_90."                         "un_Tau"                       
[127] "Tau_scf"                       "un_H"                         
[129] "H_scf"                         "un_LE"                        
[131] "LE_scf"                        "un_co2_flux"                  
[133] "co2_scf"                       "un_h2o_flux"                  
[135] "h2o_scf"                       "spikes_hf"                    
[137] "amplitude_resolution_hf"       "drop_out_hf"                  
[139] "absolute_limits_hf"            "skewness_kurtosis_hf"         
[141] "skewness_kurtosis_sf"          "discontinuities_hf"           
[143] "discontinuities_sf"            "timelag_hf"                   
[145] "timelag_sf"                    "attack_angle_hf"              
[147] "non_steady_wind_hf"            "u_spikes"                     
[149] "v_spikes"                      "w_spikes"                     
[151] "ts_spikes"                     "co2_spikes"                   
[153] "h2o_spikes"                    "chopper_LI.7500"              
[155] "detector_LI.7500"              "pll_LI.7500"                  
[157] "sync_LI.7500"                  "mean_value_RSSI_LI.7500"      
[159] "u_var"                         "v_var"                        
[161] "w_var"                         "ts_var"                       
[163] "co2_var"                       "h2o_var"                      
[165] "w.ts_cov"                      "w.co2_cov"                    
[167] "w.h2o_cov"                     "vin_sf_mean"                  
[169] "co2_mean"                      "h2o_mean"                     
[171] "dew_point_mean"                "co2_signal_strength_7500_mean"

#Add back time stamp to the combo_master_ed_met#

combo_master_ed_met$TIMESTAMP<-strptime(paste(combo_master_ed_met$date,combo_master_ed_met$time,sep=" "),format="%m/%d/%Y %H:%M", tz = "GMT")

combo_master_ed_met$TIMESTAMP[1:10]

 [1] "2019-06-25 09:00:00 GMT" "2019-06-25 09:30:00 GMT" "2019-06-25 10:00:00 GMT"
 [4] "2019-06-25 10:30:00 GMT" "2019-06-25 11:00:00 GMT" "2019-06-25 11:30:00 GMT"
 [7] "2019-06-25 12:00:00 GMT" "2019-06-25 12:30:00 GMT" "2019-06-25 13:00:00 GMT"
[10] "2019-06-25 13:30:00 GMT"

combo_master_ed_met

NA
NA
NA

#Creating a CSV File of my combine Master File!#

write.csv(combo_master_ed_met,
          paste(path.out_eddy,ver,"combo_master_ed_met",sep=""),
          quote = T,
          row.names = F)

#filtering data for quality control filters# filtering latent heat, sensible heat, co2 flux, qc_tau and h20flux by quality controls

combo_master_ed_met$LE.[!is.na(combo_master_ed_met$qc_LE)&combo_master_ed_met$qc_LE==2]<-NA

combo_master_ed_met$qqc_h2o_flux[!is.na(combo_master_ed_met$qc_h2o_flux)&combo_master_ed_met$qc_h2o_flux==2]<-NA

combo_master_ed_met$H[!is.na(combo_master_ed_met$qc_H)&combo_master_ed_met$qc_H==2]<-NA

combo_master_ed_met$u.[!is.na(combo_master_ed_met$qc_Tau)&combo_master_ed_met$qc_Tau==2]<-NA

combo_master_ed_met$co2_flux[!is.na(combo_master_ed_met$co2_flux)&combo_master_ed_met$co2_flux==2]<-NA

latent heat, sensible heat, h20_flux, and relative humidity by day for whole period

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$LE)

hist(combo_master_ed_met$LE)

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$h2o_flux)

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$H)

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$RH)

NA
NA
NA
NA
NA

#Half Hour Flux Averages

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$co2_flux,xlab='Time', ylab='Co2 Fluxes', main=' Half Hour Co2 Flux Averages',  ylim = c(-20,20),pch=1,col="blue",cex=0.4)

#Daily average of fluxes

plot(tapply(combo_master_ed_met$co2_flux ,round(combo_master_ed_met$DOY),function(x) mean(x,na.rm=T)), xlab='Time', ylab='Co2 Fluxes', main='Daily Co2 Flux Averages',pch=1,col="red",cex=1.5)

#Plotting Co2 flux by timstamp on x-axis. First half-hourly data and daily data stacked

layout(matrix(c(1,1,2,2), 2, 2, byrow = TRUE))
plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$co2_flux, xlab='Time', ylab='Co2 Fluxes', main='Half Hour Flux Averages',  ylim = c(-15,15),pch=1,col="blue",cex=0.4)
plot(tapply(combo_master_ed_met$co2_flux ,round(combo_master_ed_met$DOY),function(x) mean(x,na.rm=T)), xlab='Time', ylab='Co2 Fluxes', main='Daily Co2 Flux Averages',pch=1,col="red",cex=1.5)

NA
NA

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$co2_flux,
  ylim = c(-10, 10),
               cex=0.6,col="grey",
               xlab= 'Time', 
               ylab="", 
               main='',las=1,
               lpars=list(lwd=3,col="red"))
mtext(side=2,expression(CO[2]~Flux~'('~mu~mol~m^{-2}~s^{-1}~')'),line=2.5)

par(new= TRUE)
plot(tapply(combo_master_ed_met$co2_flux ,round(combo_master_ed_met$DOY),function(x) mean(x,na.rm=T)),
     xaxt='n', 
     xlab='Time', 
     ylab='', 
     main='', las =1,
     ylim = c(-10,10),
     lty=1,col="red",
     lwd=2,type="l")

mtext(side=2,expression(CO[2]~Flux~'('~mu~mol~m^{-2}~s^{-1}~')'),line=2.5)
abline(h=0,col="black")

#Plotting daily flux average over half hour fluxes

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$co2_flux, 
     xlab='Time', 
     ylab='Co2 Fluxes', 
     main='',  ylim = c(-10,10),pch=1,col="grey",cex=0.4)
par(new= TRUE)

plot(tapply(combo_master_ed_met$co2_flux ,round(combo_master_ed_met$DOY),function(x) mean(x,na.rm=T)),xaxt='n', xlab='Time', ylab='Co2 Fluxes', main='',ylim = c(-10,10),lty=1,col="red",lwd=2,type="l")
abline(h=0,col="black")

legend(x='bottomright',legend=c('1/2 hour fluxes', 'daily flux averages'),
col=c('grey', 'red'), pch=c(1,19))

#Daily average of Sensible Heat

plot(tapply(combo_master_ed_met$H,round(combo_master_ed_met$DOY),function(x) mean(x,na.rm=T)),xlab='Time', ylab='Sensible Heat', main='Daily Average of Sensible Heat')

#Daily average of Latent Heat

plot(tapply(combo_master_ed_met$LE ,round(combo_master_ed_met$DOY),function(x) mean(x,na.rm=T)),xlab='Time', ylab='Latent Heat', main='Daily Average of Latent Heat')

#line plot of c02_fluxes#

ggplot(data = combo_master_ed_met) + 
  geom_line(mapping = aes(x = time.id , y = co2_flux)) 

#latent heat and Sensible heat plotted over eachother

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$LE,xlab='Time', ylab='Latent and Sensible Heat', main='Latent and Sensible Heat',pch=1,col="blue",cex=0.6)
points(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$H,col="red",pch=1,cex=0.6)

legend(x='bottomright',legend=c('Latent Heat', 'Sensible Heat'),
col=c('blue', 'red'), pch=c(1,19))

#air temperatue by sensible heat an U* by co2 flux

plot(combo_master_ed_met$air_temperature-273.15,combo_master_ed_met$H,xlim=c(10,40),ylim=c(-100,500))

plot(combo_master_ed_met$u. ,combo_master_ed_met$co2_flux, ylim = c(0, 10), xlim = c(0,1))

#Adding best fit line to air temperature by sensible heat

ggplot(data = combo_master_ed_met) + 
  geom_point(mapping = aes(x = air_temperature-273.15, y = H)) +
  geom_smooth(mapping = aes(x = air_temperature-273.15, y = H))

NA
NA

#best fit line to U* bt co2 flux

ggplot(data = combo_master_ed_met) + 
  geom_point(mapping = aes(x = u., y = co2_flux)) +
  geom_smooth(mapping = aes(x = u., y = co2_flux))

library(tidyverse)
hist(combo_master_ed_met$wind_dir, xlim = c(0,370), breaks = 36, main = "Wind Direction at Concord Tower", xlab = 'Degrees')

#Comparing MET Temperature data to temperature reading from Licor intruments#

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$PTemp_C_Avg)
points(combo_master_ed_met$TIMESTAMP,combo_master_ed_met$air_temperature-273.15,col="red",pch=2)

#Adding coefficient to Net radiation#

combo_master_ed_met$Correct_NR = (combo_master_ed_met$NR_Wm2_Avg*10)/14.2

#plotting air temperature from data logger. see where our gaps line up. lines up with net radiation and gaps

plot(combo_master_ed_met$AirT_Avg)

#plotting Net Radation to see the bad data.

summary(combo_master_ed_met$Correct_NR[c(1:600,2500:4000)])

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
-106.27  -69.79  -17.23  121.53  339.30  644.37 

plot(combo_master_ed_met$Correct_NR)

#filtering out bad NR numbers

combo_master_ed_met$Correct_NR[!is.na(combo_master_ed_met$Correct_NR)&(combo_master_ed_met$Correct_NR<(-150)|combo_master_ed_met$Correct_NR>800)]<-NA
summary(combo_master_ed_met$Correct_NR)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
-146.97  -67.96  -49.57   88.05  250.14  775.35    3106 

plot(combo_master_ed_met$Correct_NR)

summary(combo_master_ed_met$Correct_NR)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
-146.97  -67.96  -49.57   88.05  250.14  775.35    3106 

#Adding coefficient to Soil Heat Flux Data

plot(combo_master_ed_met$SHF_1_mV_Avg,ylim=c(-10,10))

combo_master_ed_met$Correct_shf = (combo_master_ed_met$SHF_1_mV_Avg*16.455)
plot(combo_master_ed_met$Correct_shf,ylim=c(-50,50))

#plotting Soil heat flux to see the bad data.

summary(combo_master_ed_met$Correct_shf)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
-18.759  -9.462  -4.690  -0.672   6.417  47.868    3814 

summary(combo_master_ed_met$Correct_shf[c(2500:4000)])

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
-14.349  -8.919  -3.538   1.246  10.877  33.206 

plot(combo_master_ed_met$Correct_shf)

plot(combo_master_ed_met$Correct_shf[c(2500:4000)], ylim = c(-50,50),xlab='8/14/2019 to 9/20/2019', ylab='Soil Heat Flux ', main='Half Hour Averages of Soil Heat Flux of the Concord Site')

#filtering out bad soil heat flux numbers

combo_master_ed_met$Correct_shf[!is.na(combo_master_ed_met$Correct_shf)&(combo_master_ed_met$Correct_shf<(-20)|combo_master_ed_met$Correct_shf>50)]<-NA
summary(combo_master_ed_met$Correct_shf)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
-18.759  -9.462  -4.690  -0.672   6.417  47.868    3814 

plot(combo_master_ed_met$Correct_shf)

Energy Balance of non-gapfilled data. Slope of line is energy balance closure. Ideally it should be 1:1 Net radiation - soil heat flux= to Latent heat +sensible heat.

combo_master_ed_met$E_ng = (combo_master_ed_met$Correct_NR-combo_master_ed_met$Correct_shf)

summary(combo_master_ed_met$E_ng )

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
-138.33  -60.82  -43.80   75.92  206.58  763.13    3901 

plot(combo_master_ed_met$E_ng)

combo_master_ed_met$E_le_and_H =(combo_master_ed_met$LE+combo_master_ed_met$H)
summary(combo_master_ed_met$E_le_and_H  )

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
-514.86  -18.34    3.21   74.61  164.30  523.76     882 

plot(combo_master_ed_met$E_le_and_H )

scatter.smooth(combo_master_ed_met$E_ng, combo_master_ed_met$E_le_and_H )

lm(combo_master_ed_met$E_le_and_H  ~ combo_master_ed_met$E_ng)

Dropping 4195 rows with missing values

Call:
lm(formula = combo_master_ed_met$E_le_and_H ~ combo_master_ed_met$E_ng)

Coefficients:
             (Intercept)  combo_master_ed_met$E_ng  
                 16.9654                    0.5162  

summary(lm(combo_master_ed_met$E_le_and_H  ~ combo_master_ed_met$E_ng-1))

Dropping 4195 rows with missing values

Call:
lm(formula = combo_master_ed_met$E_le_and_H ~ combo_master_ed_met$E_ng - 
    1)

Residuals:
    Min      1Q  Median      3Q     Max 
-435.88    3.73   16.47   31.23  285.85 

Coefficients:
                         Estimate Std. Error t value Pr(>|t|)    
combo_master_ed_met$E_ng 0.546294   0.004323   126.4   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 50.35 on 2817 degrees of freedom
  (4195 observations deleted due to missingness)
Multiple R-squared:   0.85, Adjusted R-squared:   0.85 
F-statistic: 1.597e+04 on 1 and 2817 DF,  p-value: < 2.2e-16

#1/2 hour and daily averages of Latent Heat

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$LE, 
     xlab='Time', 
     ylab= expression(Latent~Heat~'('~W~m^{-2}~')'), 
     main='',  
     ylim = c(-40,300),
     pch=1,
     col="grey",cex=0.4)
par(new= TRUE)

plot(tapply(combo_master_ed_met$LE ,round(combo_master_ed_met$DOY),function(x) mean(x,na.rm=T)),
     xaxt='n', 
     xlab='Time', 
    ylab= expression(Latent~Heat~'('~W~m^{-2}~')'),  
     main='',
     ylim = c(-40,300),lty=1,col="red",lwd=2,type="l")
abline(h=0,col="black")

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$LE,
  ylim = c(-40, 300),
               cex=0.4,col="grey",
               xlab= 'Time', 
               ylab="", 
               main='',las=1,
               lpars=list(lwd=3,col="red"))
mtext(side=2,expression(Latent~Heat~'('~W~m^{-2}~')'),line=2.5)

par(new= TRUE)
plot(tapply(combo_master_ed_met$LE ,round(combo_master_ed_met$DOY),function(x) mean(x,na.rm=T)),
     xaxt='n', 
     xlab='Time', 
     ylab='', 
     main='', las =1,
     ylim = c(-40,300),
     col="red",pch=1,cex=1.0)

mtext(side=2,expression(Latent~Heat~'('~W~m^{-2}~')'),line=2.5)
abline(h=0,col="black")

legend(x='topright',legend=c('1/2 hour averages', 'daily averages'),
col=c('grey', 'red'), pch=c(19,1))

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$LE, 
     xlab='Time', 
     ylab='Latent Heat', 
     main='Latent Heat: Daily and Half Hour Averages', 
     ylim = c(-40,300),pch=1,col="grey",cex=0.4)
par(new= TRUE)

plot(tapply(combo_master_ed_met$LE ,round(combo_master_ed_met$DOY),function(x) mean(x,na.rm=T)),xaxt='n'
     , xlab='Time',
     ylab='Latent Heat', 
     main='Latent Heat: Daily and Half Hour Averages',
     ylim = c(-40,300),
     col="red",pch=1,cex=1.0)
abline(h=0,col="black")

legend(x='topright',legend=c('1/2 hour averages', 'daily averages'),
col=c('grey', 'red'), pch=c(1,1))

#1/2 hour and daily averages of Sensible Heat

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$H, xlab='Time', ylab='Sensible Heat', main='Sensible Heat: Daily and Half Hour Averages',  ylim = c(-40,300),pch=1,col="grey",cex=0.4)
par(new= TRUE)

plot(tapply(combo_master_ed_met$H ,round(combo_master_ed_met$DOY),function(x) mean(x,na.rm=T)),xaxt='n', xlab='Time', ylab='Sensible Heat', main='Sensible Heat: Daily and Half Hour Averages',ylim = c(-40,300),lty=1,col="red",lwd=2,type="l")
abline(h=0,col="black")

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$H, xlab='Time', ylab='Sensible Heat', main='Sensible Heat: Daily and Half Hour Averages',  ylim = c(-40,300),pch=1,col="grey",cex=0.4)
par(new= TRUE)

plot(tapply(combo_master_ed_met$H ,round(combo_master_ed_met$DOY),function(x) mean(x,na.rm=T)),xaxt='n', xlab='Time', ylab='Sensible Heat', main='Sensible Heat: Daily and Half Hour Averages',ylim = c(-40,300),col="red",pch=1,cex=1.0)
abline(h=0,col="black")

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$H,
  ylim = c(-40, 300),
               cex=0.4,col="grey",
               xlab= 'Time', 
               ylab="", 
               main='',las=1,
               lpars=list(lwd=3,col="red"))
mtext(side=2,expression(Sensible~Heat~'('~W~m^{-2}~')'),line=2.5)

par(new= TRUE)
plot(tapply(combo_master_ed_met$H ,round(combo_master_ed_met$DOY),function(x) mean(x,na.rm=T)),
     xaxt='n', 
     xlab='Time', 
     ylab='', 
     main='', las =1,
     ylim = c(-40,300),
     col="red",pch=1,cex=1.0)

mtext(side=2,expression(Sensible~Heat~'('~W~m^{-2}~')'),line=2.5)
abline(h=0,col="black")

#Net Radiation Daily and 1/2 hour averages

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$Correct_NR, xlab='Time', ylab='Net Radiation', main='Net Radiation: Daily and Half Hour Averages',  ylim = c(-80,800),pch=1,col="grey",cex=0.4)
par(new= TRUE)

plot(tapply(combo_master_ed_met$Correct_NR ,round(combo_master_ed_met$DOY),function(x) mean(x,na.rm=T)),xaxt='n', xlab='Time', ylab='Net Radiation', main='Net Radiation: Daily and Half Hour Averages',ylim = c(-80,800),lty=1,col="red",lwd=2,type="l")
abline(h=0,col="black")

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$Correct_NR, xlab='Time', ylab='Net Radiation', main='Net Radiation: Daily and Half Hour Averages',  ylim = c(-80,800),pch=1,col="grey",cex=0.4)
par(new= TRUE)

plot(tapply(combo_master_ed_met$Correct_NR ,round(combo_master_ed_met$DOY),function(x) mean(x,na.rm=T)),xaxt='n', xlab='Time', ylab='Net Radiation', main='Net Radiation: Daily and Half Hour Averages',ylim = c(-80,800),col="red",pch=1,cex=1.0)
abline(h=0,col="black")

#soil heat flux

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$Correct_shf, xlab='Time', ylab='Soil Heat Flux', main='Soil Heat Flux: Daily and Half Hour Averages',  ylim = c(-20,40),pch=1,col="grey",cex=0.4)
par(new= TRUE)

plot(tapply(combo_master_ed_met$Correct_shf,round(combo_master_ed_met$DOY),function(x) mean(x,na.rm=T)),xaxt='n', xlab='Time', ylab='Soil Heat Flux', main='Soil Heat Flux: Daily and Half Hour Averages',ylim = c(-20,40),lty=1,col="red",lwd=2,type="l")
abline(h=0,col="black")

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$Correct_shf, xlab='Time', ylab='Soil Heat Flux', main='Soil Heat Flux: Daily and Half Hour Averages',  ylim = c(-20,40),pch=1,col="grey",cex=0.4)
par(new= TRUE)

plot(tapply(combo_master_ed_met$Correct_shf,round(combo_master_ed_met$DOY),function(x) mean(x,na.rm=T)),xaxt='n', xlab='Time', ylab='Soil Heat Flux', main='Soil Heat Flux: Daily and Half Hour Averages',ylim = c(-20,40),col="red",pch=1,cex=1.0)
abline(h=0,col="black")

#air temperature daily and 1/2 hour averages

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$air_temperature-273.15, xlab='Time', ylab='Air Temperature C: Daily and 1/2 Hour Averages', main=' ',  ylim = c(5,40),pch=1,col="grey",cex=0.4)
par(new= TRUE)

plot(tapply(combo_master_ed_met$air_temperature-273.15,round(combo_master_ed_met$DOY),function(x) mean(x,na.rm=T)),xaxt='n', xlab='Time', ylab='Air Temperature C: Daily and 1/2 Hour Averages', main='',ylim = c(5,40),lty=1,col="red",lwd=2,type="l")
abline(h=0,col="black")

#Water Fluxes daily and 1/2 hour averages

plot(combo_master_ed_met$TIMESTAMP ,combo_master_ed_met$h2o_flux, xlab='Time', ylab='H2O Fluxes', main='H2O Fluxes: Daily and Half Hour Averages',  ylim = c(-2,5),pch=1,col="grey",cex=0.4)
par(new= TRUE)

plot(tapply(combo_master_ed_met$h2o_flux,
            round(combo_master_ed_met$DOY),
            function(x) mean(x,na.rm=T)),
     xaxt='n',
     xlab='Time',
     ylab='H2O Fluxes',
     main='H2O Fluxes: Daily and Half Hour Averages',
     ylim = c(-2,5),
     lty=1,
     col="red",
     lwd=2,
     type="l")
abline(h=0,col="black")

plot(cumsum(tapply(combo_master_ed_met$h2o_flux,
            round(combo_master_ed_met$DOY),
            function(x) mean(x,na.rm=T)))*18.02/1000000*1800*48,
     xaxt='n',
     xlab='Days since June 25',
     ylab=expression(Cumulative~Evapotranspiration~'('~mm~')'),
     main='',
     ylim = c(0,150),
     lty=1,
     col="red",
     lwd=2,
     type="l")
axis(side=1,at=seq(0,120,by=30))

plot(cumsum(tapply(combo_master_ed_met$co2_flux,
            round(combo_master_ed_met$DOY),
            function(x) mean(x,na.rm=T)))*12/1000000*1800*48,
     xaxt='n',
     xlab='Days since June 25',
     ylab=expression(Cumulative~NEE~'('~g~C~m^{-2}~')'),
     main='',
     ylim = c(-50,10),
     lty=1,
     col="red",
     lwd=2,
     type="l")
axis(side=1,at=seq(0,120,by=30))

#Co2 fluxes v.s u star, temperature, and VPD

plot(combo_master_ed_met$u. ,combo_master_ed_met$co2_flux, ylim = c(-5, 10), xlim = c(0,1), xlab='U*', ylab='Co2 Fluxes', main='Co2 Fluxes v.s U*')

plot(combo_master_ed_met$VPD ,combo_master_ed_met$co2_flux, ylim = c(-5, 10), xlab='VPD', ylab='Co2 Fluxes', main='Co2 Fluxes v.s VPD*')

plot(combo_master_ed_met$air_temperature-273.15 ,combo_master_ed_met$co2_flux, ylim = c(-5, 10), xlab='Air Temperature (C)', ylab='Co2 Fluxes', main='Co2 Fluxes v.s Air temperature*')

plot(combo_master_ed_met$Correct_NR[combo_master_ed_met$daytime==1] ,combo_master_ed_met$co2_flux[combo_master_ed_met$daytime==1], ylim = c(-10, 10), xlab='Net Radiation', ylab='Co2 Fluxes', main='Co2 Fluxes v.s Net Radiation')

summary(combo_master_ed_met$VPD)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
    0.0   436.5   947.7  1236.6  1785.0  5784.0     317 

scatter.smooth(combo_master_ed_met$Correct_NR[combo_master_ed_met$daytime==1],
               combo_master_ed_met$co2_flux[combo_master_ed_met$daytime==1],
               ylim = c(-10, 10),
               cex=0.6,col="grey",
               xlab=expression(Daytime~Net~Radiation~'('~W~m^{-2}~')'), 
               ylab="", 
               main='',las=1,
               lpars=list(lwd=3,col="red"))
mtext(side=2,expression(CO[2]~Flux~'('~mu~mol~m^{-2}~s^{-1}~')'),line=2.5)

scatter.smooth(combo_master_ed_met$air_temperature[combo_master_ed_met$daytime==0&
                                                     combo_master_ed_met$u.>=0.15]-273.15,
               combo_master_ed_met$co2_flux[combo_master_ed_met$daytime==0&
                                                     combo_master_ed_met$u.>=0.15],
               ylim = c(-5, 10),
               cex=0.6,col="grey",
               xlab=expression(Air~Temperature~'('~degree~C~')'), 
               ylab="", 
               main='',las=1,
               lpars=list(lwd=3,col="red"))
mtext(side=2,expression(CO[2]~Flux~'('~mu~mol~m^{-2}~s^{-1}~')'),line=2.5)

#cumulative